/*  
 * Copyright LWJGL. All rights reserved.    
 * License terms: https://www.lwjgl.org/license 
 */
package org.lwjgl.demo.util;

import java.util.ArrayList;
import java.util.List;
import org.joml.Vector3i;

/**
 * Custom marching cubes algorithm that generates vertices on an integer lattice
 * without linear interpolating the density values.
 * 
 * @author Kai Burjack
 */
public class MarchingCubes {

    private static final short[] E = { 0x0, 0x109, 0x203, 0x30a, 0x406, 0x50f, 0x605, 0x70c, 0x80c, 0x905, 0xa0f, 0xb06,
            0xc0a, 0xd03, 0xe09, 0xf00, 0x190, 0x99, 0x393, 0x29a, 0x596, 0x49f, 0x795, 0x69c, 0x99c, 0x895, 0xb9f,
            0xa96, 0xd9a, 0xc93, 0xf99, 0xe90, 0x230, 0x339, 0x33, 0x13a, 0x636, 0x73f, 0x435, 0x53c, 0xa3c, 0xb35,
            0x83f, 0x936, 0xe3a, 0xf33, 0xc39, 0xd30, 0x3a0, 0x2a9, 0x1a3, 0xaa, 0x7a6, 0x6af, 0x5a5, 0x4ac, 0xbac,
            0xaa5, 0x9af, 0x8a6, 0xfaa, 0xea3, 0xda9, 0xca0, 0x460, 0x569, 0x663, 0x76a, 0x66, 0x16f, 0x265, 0x36c,
            0xc6c, 0xd65, 0xe6f, 0xf66, 0x86a, 0x963, 0xa69, 0xb60, 0x5f0, 0x4f9, 0x7f3, 0x6fa, 0x1f6, 0xff, 0x3f5,
            0x2fc, 0xdfc, 0xcf5, 0xfff, 0xef6, 0x9fa, 0x8f3, 0xbf9, 0xaf0, 0x650, 0x759, 0x453, 0x55a, 0x256, 0x35f,
            0x55, 0x15c, 0xe5c, 0xf55, 0xc5f, 0xd56, 0xa5a, 0xb53, 0x859, 0x950, 0x7c0, 0x6c9, 0x5c3, 0x4ca, 0x3c6,
            0x2cf, 0x1c5, 0xcc, 0xfcc, 0xec5, 0xdcf, 0xcc6, 0xbca, 0xac3, 0x9c9, 0x8c0, 0x8c0, 0x9c9, 0xac3, 0xbca,
            0xcc6, 0xdcf, 0xec5, 0xfcc, 0xcc, 0x1c5, 0x2cf, 0x3c6, 0x4ca, 0x5c3, 0x6c9, 0x7c0, 0x950, 0x859, 0xb53,
            0xa5a, 0xd56, 0xc5f, 0xf55, 0xe5c, 0x15c, 0x55, 0x35f, 0x256, 0x55a, 0x453, 0x759, 0x650, 0xaf0, 0xbf9,
            0x8f3, 0x9fa, 0xef6, 0xfff, 0xcf5, 0xdfc, 0x2fc, 0x3f5, 0xff, 0x1f6, 0x6fa, 0x7f3, 0x4f9, 0x5f0, 0xb60,
            0xa69, 0x963, 0x86a, 0xf66, 0xe6f, 0xd65, 0xc6c, 0x36c, 0x265, 0x16f, 0x66, 0x76a, 0x663, 0x569, 0x460,
            0xca0, 0xda9, 0xea3, 0xfaa, 0x8a6, 0x9af, 0xaa5, 0xbac, 0x4ac, 0x5a5, 0x6af, 0x7a6, 0xaa, 0x1a3, 0x2a9,
            0x3a0, 0xd30, 0xc39, 0xf33, 0xe3a, 0x936, 0x83f, 0xb35, 0xa3c, 0x53c, 0x435, 0x73f, 0x636, 0x13a, 0x33,
            0x339, 0x230, 0xe90, 0xf99, 0xc93, 0xd9a, 0xa96, 0xb9f, 0x895, 0x99c, 0x69c, 0x795, 0x49f, 0x596, 0x29a,
            0x393, 0x99, 0x190, 0xf00, 0xe09, 0xd03, 0xc0a, 0xb06, 0xa0f, 0x905, 0x80c, 0x70c, 0x605, 0x50f, 0x406,
            0x30a, 0x203, 0x109, 0x0 };

    private static final long[] TS = { 0x0L, 0x491L, 0xa21L, 0x29a492L, 0xb32L, 0xb32491L, 0xa31b3aL, 0x9ab9b3493L,
            0x3c4L, 0x1c93c1L, 0xc431a2L, 0xc9aca23c2L, 0x4bc2b4L, 0xbc9b912b1L, 0xabcac41a4L, 0xc9bb9aL, 0x985L,
            0x548145L, 0x859a21L, 0x248285a25L, 0x859b32L, 0xb32514854L, 0x85931ab3aL, 0x5a84838a3ab3L, 0x3c4859L,
            0x51353c85cL, 0xc4385921aL, 0x23a3cac5ac85L, 0x598bc42b4L, 0x5c8512c52bc2L, 0x41cbcac1a985L, 0xbcaac5c85L,
            0x56aL, 0x49156aL, 0x162561L, 0x624649569L, 0x56ab32L, 0x6a5b32914L, 0x315356b36L, 0x9545646346b3L,
            0xc4356aL, 0x6a5c913c1L, 0xc43621561L, 0x695c93963623L, 0x56a42bc4bL, 0xbc92b92916a5L, 0x41cbc6c16156L,
            0xc9bb96956L, 0xa8698aL, 0x48646a14aL, 0x862821981L, 0x864462L, 0x32b86a98aL, 0x48614616a32bL,
            0x36b869639319L, 0x8644636b3L, 0x3c4a986a8L, 0xc8313a38a86aL, 0x862982921c43L, 0x62882c23cL,
            0xc4b42b86996aL, 0x1bcb121c8a16186L, 0x18681916b41c1bcL, 0x6c86bcL, 0x67bL, 0x7b6491L, 0x7b621aL,
            0x7b69a2492L, 0x273672L, 0x914732672L, 0x73171a67aL, 0x9347363969a6L, 0x67bc43L, 0x67b13c91cL, 0x7b6c43a21L,
            0xc9a3ca3a27b6L, 0x426467c47L, 0x7c62616c1c91L, 0xa616717417c4L, 0x9acca7a67L, 0x9857b6L, 0xb67485145L,
            0x8597b61a2L, 0x5a84828a267bL, 0x985267327L, 0x854514736632L, 0x73167161a859L, 0xa737a6a345a8a48L,
            0x67b5983c4L, 0xc831353857b6L, 0x7b6c43985a21L, 0x7b65c8c5a3ca23aL, 0x7c62646c4859L, 0xc515c8c127c6c26L,
            0x85947c741671a61L, 0xac8a85ca7a67L, 0xb57a5bL, 0x491ba57b5L, 0x15717b21bL, 0xb27579729249L, 0x573532a52L,
            0x573a53a32914L, 0x735531L, 0x735539349L, 0x43c57ba5bL, 0x7b5ba5c93391L, 0xb275717213c4L,
            0x2c9c23295b27257L, 0x47c42a74a57aL, 0x25752a27c1292c9L, 0x5711747c4L, 0x97c957L, 0xba9b987b8L,
            0xb87ba18b1481L, 0x1929828b287bL, 0x48228b87bL, 0x879a92972732L, 0xa484a1a872a3a73L, 0x317718198L,
            0x387348L, 0x879a9b97bc43L, 0x8bab878a1c83813L, 0xc43b878b2982192L, 0x28727b82c23cL, 0x74247c72a8797a9L,
            0x87c2a1L, 0x17c1c4718198L, 0x7c8L, 0xc78L, 0x78c914L, 0x78ca21L, 0x78c249a29L, 0x8c732bL, 0x8c7914b32L,
            0x8c7ab31a3L, 0x9ab49b4b38c7L, 0x837438L, 0x137178918L, 0xa21843783L, 0x7899a2792372L, 0x84282b78bL,
            0x912892b8278bL, 0x8b7ab1b81841L, 0xab99b8b78L, 0x79c597L, 0x751714c74L, 0x21a759c79L, 0x74c24a47a75aL,
            0x2b39c7597L, 0x751c71c14b32L, 0xab3a31c759c5L, 0x47574c45a34b4abL, 0x375359439L, 0x375351L,
            0x9457535431a2L, 0x7533525a2L, 0x2b7759279429L, 0x51771b12bL, 0x4aba414b7945475L, 0x5b75abL, 0xc786a5L,
            0x78c6a5491L, 0xc78156216L, 0x62456454978cL, 0xc7832b56aL, 0x6a5491b328c7L, 0x315b35b56c78L,
            0x78c36b634564954L, 0xa56378438L, 0x89737179156aL, 0x156162784374L, 0x962695923897937L, 0x84278272b56aL,
            0x56a189812782b72L, 0xb848b7b416b5b15L, 0xb95b569b8b78L, 0xa9cac76a7L, 0x6a1761471c74L, 0xc76621c619c1L,
            0x2466474c7L, 0x76c9cac6ab32L, 0xb32a767a1c714c1L, 0x63136b61976c69cL, 0x46b4b36474c7L, 0x394376936a96L,
            0x37117a76aL, 0x937394976192962L, 0x723762L, 0x7a9a76794b72742L, 0x17616a71b12bL, 0xb76941L, 0x76bL,
            0xc68b6cL, 0x14968cb6cL, 0x1a2cb68c6L, 0x24929a8cb68bL, 0x26828c32cL, 0xc38682832491L, 0x8c331a83a68aL,
            0x39a9343a6c38368L, 0x684643b63L, 0x63b689369139L, 0x3b48464b621aL, 0x36863b38923a39aL, 0x684642L,
            0x682281891L, 0x84664a41aL, 0x8a689aL, 0x9cb9b6596L, 0x14ccb61c6516L, 0x65bcb9b59a21L, 0x52425a54c65b5cbL,
            0x9659c3693263L, 0xc262c3c654c1c51L, 0x69c9656c3a61631L, 0x4c365aL, 0x594654364b63L, 0x1355363b6L,
            0xa21965694b643b4L, 0x35a3a25363b6L, 0x264469659L, 0x612651L, 0x64161a469659L, 0x65aL, 0xcbaca58c5L,
            0xcba8ca8a5491L, 0xc581525c2cb2L, 0x5cbc585b2954524L, 0x32ac3a5ca8c5L, 0x4912c3c2a8ca58aL, 0x15335c58cL,
            0x53454935c58cL, 0xa58843a83ba3L, 0x81318983b58a8baL, 0xb151b2b583b4b84L, 0x5893b2L, 0x4288252a5L,
            0x2892918252a5L, 0x458415L, 0x895L, 0x9cb9baL, 0xbacca4a14L, 0xcb99b1b21L, 0xb4cb24L, 0x9caac2c32L,
            0xac3a32ca4a14L, 0xc19c31L, 0xc34L, 0xa9bb93943L, 0x3a13baL, 0x9b2921b93943L, 0x3b2L, 0x92a942L, 0x2a1L,
            0x941L, 0x0L };

    private static final byte[] CS = { 0x1, 0x6, 0x9, 0x4, 0x21, 0x26, 0x29, 0x24, 0x10, 0x12, 0x1a, 0x18 };

    public static List<Vector3i> march(byte[] ds, byte iso, int dx, int dy, int dz) {
        List<Vector3i> vts = new ArrayList<>();
        int vl[] = new int[12 * 3];
        for (int z = 0; z < dz - 1; z++) {
            for (int y = 0; y < dy - 1; y++) {
                for (int x = 0; x < dx - 1; x++) {
                    int ci = determineCase(ds, iso, dx, dy, z, y, x);
                    int b = E[ci];
                    if (b == 0)
                        continue;
                    generateVertices(vl, z, y, x, b);
                    addVertices(vts, vl, TS[ci]);
                }
            }
        }
        return vts;
    }

    private static void addVertices(List<Vector3i> vts, int[] vl, long c) {
        while (c != 0L) {
            int i1 = (int) (c & 0xF) - 1, i2 = (int) (c >>> 4 & 0xF) - 1, i3 = (int) (c >>> 8 & 0xF) - 1;
            vts.add(new Vector3i(vl[3 * i3], vl[3 * i3 + 1], vl[3 * i3 + 2]));
            vts.add(new Vector3i(vl[3 * i2], vl[3 * i2 + 1], vl[3 * i2 + 2]));
            vts.add(new Vector3i(vl[3 * i1], vl[3 * i1 + 1], vl[3 * i1 + 2]));
            c >>>= 12;
        }
    }

    private static void generateVertices(int[] vl, int z, int y, int x, int b) {
        for (int i = 0; i < 12; i++) {
            if ((b & 1 << i) == 0)
                continue;
            vl[3 * i] = (x << 1) + (CS[i] & 3);
            vl[3 * i + 1] = (y << 1) + (CS[i] >>> 2 & 3);
            vl[3 * i + 2] = (z << 1) + (CS[i] >>> 4 & 3);
        }
    }

    private static int determineCase(byte[] ds, byte iso, int dx, int dy, int z, int y, int x) {
        int x0y0z0 = x + dx * y + dx * dy * z, x1y0z0 = x0y0z0 + 1, x0y1z0 = x0y0z0 + dx, x1y1z0 = x0y1z0 + 1,
                x0y0z1 = x0y0z0 + dx * dy, x1y0z1 = x1y0z0 + dx * dy, x0y1z1 = x0y1z0 + dx * dy,
                x1y1z1 = x1y1z0 + dx * dy;
        int ci = 0;
        ci |= ds[x0y0z0] > iso ? 1 << 0 : 0;
        ci |= ds[x1y0z0] > iso ? 1 << 1 : 0;
        ci |= ds[x1y1z0] > iso ? 1 << 2 : 0;
        ci |= ds[x0y1z0] > iso ? 1 << 3 : 0;
        ci |= ds[x0y0z1] > iso ? 1 << 4 : 0;
        ci |= ds[x1y0z1] > iso ? 1 << 5 : 0;
        ci |= ds[x1y1z1] > iso ? 1 << 6 : 0;
        ci |= ds[x0y1z1] > iso ? 1 << 7 : 0;
        return ci;
    }

}
